HPAVC/HPAVC CD-ROM.iso/pc/SIMULATE.ZIP/simulate.pas

Pascal/Delphi Source File | 1994-08-21 | 38.0 KB | 1,273 lines

      (************************  BEGIN  **************************
      **                                                        **
      **                                                        **
      **                    Flight  Simulator                   **
      **                  Subsonic Jet Aircraft                 **
      **                      Version  4.26                     **
      **                                                        **
      **                                                        **
      ************************************************************

      ************************************************************
      *                                                          *
      *                                                          *
      *            Captain Larry Thomas Brewster, Ph.D.          *
      *                                                          *
      *                University of Missouri-Rolla              *
      *                Computer Science Department               *
      *                       Feb 26, 1990                       *
      *                                                          *
      *                                                          *
      ***********************************************************)


      (***********************************************************
      This software is made available only to individuals and only
      for educational purposes. Any and all commercial use is
      stricly prohibited.
      ***********************************************************)


      (***********************************************************
      This program is designed to run on a Toshiba T3100 portable
      computer using Borlands Pascal Version 4.0 or later compiler.
      Changes in the graphics routines may be necessary
      for any other computer.
      ***********************************************************)


      (********              KNOWN BUGS

      program causes math coprocessor stack overflows on an IBM XT

      bank angles of exactly zero, ninety, one-eighty, etc. cause the
      horizon bar to have a length of zero.

      angle of attack (alpha) is positive in inverted flight

      graphics routines are crude - may require modification if
      not 640 X 400

      this program is not a polished end product and as such no
      guarantees as to performance or accuracy is implied
      it is intended as an educational tool for experimentation
      ********)

{$N-} { turn off math chip to prevent stack overflow }

program Simulate;

  uses Crt, Graph;

  const

  {       Aircraft Specifications             Random Aircraft          }

    Clmax                  =      1.20;   { max coef of lift           }
    AR                     =      7.62;   { aspect ratio span^2/Sref   }
    Tail_Arm               =        25;   { dist from a/c a.c. to tail }
    max_alpha_tail         =      0.35;   { max deflection of elevator }
    k_elevator             =      3.50;   { elevator constant          }
    max_aileron            =      0.35;   { max aileron deflection     }
    k_aileron              =      1.00;   { just a guess to look right }
    maximum_mach           =      0.90;   { typical Mmo                }
    maximum_ias            =    500.00;   { max guage airspeed         }
    minimum_weight         =  22750.00;   { min inflight weight        }
    maximum_weight         =  45500.00;   { max inflight weight        }
    wing_loading           =     31.30;   { max weight over Sref       }
    thrust_to_weight       =      0.30;   { thrust to weight ratio     }
    maximum_altitude       =  55000.00;   { aircraft absolute ceiling  }
    minimum_rpm            =     0.063;
    maximum_rpm            =         1;
    max_engine_rpm         =       100;

  {                  End of Aircraft Specifications                    }

  {              Constants used by equations of flight                 }

    dtime                  =      1.13;   { seconds @8MHz 2.3754 @4MHz }
    text                   =         8;   { graphics text width&height }
    tropopause             =  36089.00;   { feet                       }
    sea_level_temperature  =   518.688;   { degrees Rankine            }
    tropopause_temperature =   389.988;   { degrees Rankine            }
    sea_level_density      = 0.0023769;   { slugs / cubic foot         }
    Mach_1                 = 49.021424;   { sqrt(1.4*1716.5)           }
    nautical_miles_to_feet =  6076.118;   { 5280*1.15078 NM to feet    }
    feet_to_knots          = 0.5924835;   { 3600/(5280*1.15078)        }
    knots_to_feet          = 1.6878107;   { (5280*1.15078)/3600        }
    seconds_to_hour        =   3600.00;   { seconds in an hour         }
    g                      = 32.174049;   { feet / sec squared         }
    minimum_altitude       =      0.00;   { density equation limit     }
    dev                    =     1E-20;   { minimum value used in plot }

  {             End of Constants used by equations of flight           }

  var

    w,          mach,     altitude,  attitude,
    trim,       theta,    psi,       alpha_tail,  alpha_aileron,
    att_last,   psi_last, tas_last,  n_last,      mach_last,
    theta_last,           alt_last,               rpm_last,
    last_fuelflow,        rpm,                    R : real;

    tas_lastx,  tas_lasty,  alt_lastx,  alt_lasty,
    xl, xm,  yl, ym, xx, yy, zz, rcx, rcy, tbx, tby,
    sky_x, sky_y, rpm_lastx, rpm_lasty,

  { graphics values }

    col_1, col_2, col_3, col_4,

    row_1, row_2, row_3, radii,

    ias_x, hdi_x, alt_x, rpm_x,
    ias_y, hdi_y, alt_y, rpm_y,
    ias_r, hdi_r, alt_r, rpm_r,

    tab_x, hsi_x, roc_x,
    tab_y, hsi_y, roc_y,
    tab_r, hsi_r, roc_r,

    GraphDriver, GraphMode, ErrorCode : integer;
    Xasp, Yasp   : word;
    gAR          : real;
    function_key : boolean;
    SaveExitProc : pointer;
    s            : string;
    ch           : char;


function sigma : real;     { density / sea_level_density }
  begin   { density-ratio }
    sigma := exp( - (altitude/23111) + 0.294 * sin(altitude/28860)
                                     + 0.213 * sin(altitude/86580) );
  end {sigma};


function density : real;   { slugs per cubic foot }
  begin
    density := sea_level_density * sigma;
  end {density};


function temperature : real; { degrees Rankine i.e. F + 460 }
  begin
    case ( altitude > tropopause ) of
      true  :  temperature := tropopause_temperature;
      false :  temperature := ( tropopause_temperature
                              - sea_level_temperature )
                              * ( altitude / tropopause )
                              + sea_level_temperature;
    end; { case }
  end; { temperature }


function mach_one : real;  { feet per second }
  begin
    mach_one := mach_1 * sqrt( temperature );
  end; { mach_one }


function V : real;  { feet per second }
  begin
    V := mach * mach_one;
  end; { V }


function ias : real;    { feet per second }
  begin
    ias := V * sqrt(sigma);
  end; { ias }


function calculate_mach_from_ias(ias_inp: real) : real;
  begin           { input in knots!  output in ft/sec }
    calculate_mach_from_ias  := ias_inp * knots_to_feet
                                / (sqrt(sigma)*mach_one);
  end {calculate_mach_from_ias};


function calculate_mach_from_tas(tas_inp: real) : real;
  begin           { input in knots!  output in ft/sec }
    calculate_mach_from_tas  := tas_inp * knots_to_feet / mach_one;
  end {calculate_mach_from_tas};


function m : real;  { slugs "pounds-sec^2/ft" }
  begin
    m := w / g;
  end; { m }


function q : real;
  begin
    q := density * sqr( V ) / 2.0;
  end; { q }


function n : real;  { load factor - or g-loading }
  begin
    n := ( alpha_tail + trim ) * k_elevator * q * tail_arm / m;
  end; { n }


function Cl : real;  { coefficient of lift }
  begin
    Cl := n * w / ( q * maximum_weight/wing_loading );
  end; { Cl }


function alpha : real;  { wing angle of attack }
  const Cla0  =  0.10;  { coefficient of lift at zero angle of attack }
        dClda =  3.53;  { slope of Cl versus alpha }
  begin
    alpha :=  (Cl - Cla0) / dClda;
  end; { alpha }


function gamma : real;  { flight path angle }
  begin
    gamma := attitude - alpha;
  end; { gamma }


function e : real;  { efficiency - from drag data }
  const k1 = 0.87000;
        k2 = 0.08355;
  begin
    e := k1 + k2 * sqr(mach);
  end {e};


function CDl : real;  { coefficient of drag due to lift }
  begin
    CDl := sqr( Cl ) / ( PI * AR * e );
  end; { CDl }


function Cd0 : real;  { coefficient of parasite drag due to airframe }
  begin
    Cd0   := 0.02;
  end; { Cd0 }


function D : real;    { total a/c drag }
  begin
    D := (Cd0 + CDl) * q * maximum_weight/wing_loading;
  end; { D }


function Xdot : real;  { horizontal velocity }
  begin
    Xdot := abs(V * cos(gamma));
  end; { Xdot }


function thetadot : real;  { bank angle rate of change }
  begin
    thetadot := alpha_aileron * k_aileron * V
                / sqrt(AR * maximum_weight/wing_loading);
  end; { thetadot }


function PsiDot : real; { horizontal turning rate }
  const  k = 0.17; { cos(80) }
  begin
    if abs( cos(gamma) ) > k then
         PsiDot := (g * n * sin(theta)) / V * cos(gamma)
    else Psidot := thetadot;
  end; { PsiDot }


function gammadot : real;  { rate of change of flight path }
  begin
    gammadot := ( n * cos(theta) - cos(gamma) ) * g / V;
  end; { gammadot }


function maximum_thrust : real;   { in pounds }
  begin
    maximum_thrust := maximum_weight * thrust_to_weight * sqrt(sigma);
  end {maximum_thrust};


function T : real;   { thrust in pounds }
  begin
    T := maximum_thrust * rpm;
  end {T};


function Edot : real;  { energy rate of change - power }
  begin
    Edot := V * (T * cos(alpha) - D);
  end; { Edot }


function Vdot : real;  { acceleration }
  begin
    Vdot := ( (T * cos(alpha) - D) / m ) - g * sin(gamma);
  end; { Vdot }


function Hdot : real;  { vertical rate of climb }
  begin
    Hdot := (Edot / w) - (V * Vdot / g);
  end; { Hdot }


function Vstall : real;  {  power off flight }
  begin
    Vstall := sqrt( abs( (w + abs(n) * w)
          / ( sea_level_density * Clmax
              * maximum_weight/wing_loading)) );
  end; { Vstall }


function fuelflow : real;       { in pounds per second           }
  begin                         { sfc = 0.25 lb/hr per lb thrust }
      fuelflow := T * 0.25 / seconds_to_hour;
  end {fuelflow};


procedure read_character;
  begin
    ch    := readkey;
    if ch <> #0 then function_key := false
    else
      begin
        function_key :=    true;
        ch           := readkey;
      end;
  end { read_character };


procedure filter_weight;
  begin
    if      w < minimum_weight then w := minimum_weight
    else if w > maximum_weight then w := maximum_weight ;
  end {filter_weight};


procedure filter_altitude;
 begin
  if altitude > maximum_altitude then altitude := maximum_altitude;
  if altitude < minimum_altitude then altitude := minimum_altitude;
 end {filter_altitude};


procedure filter_rpm;
  begin
    if      rpm  <  minimum_rpm / sqrt(sigma) then
            rpm :=  minimum_rpm / sqrt(sigma)
    else if rpm  >  maximum_rpm  then
            rpm :=  maximum_rpm;
  end {filter_rpm};


procedure filter_mach;
  begin
    if  mach > maximum_mach then mach := calculate_mach_from_ias(mach);
    if ias * feet_to_knots  >  maximum_ias then
        mach := calculate_mach_from_ias(maximum_ias)
    else if ias < Vstall then
        mach := calculate_mach_from_ias(Vstall * feet_to_knots);
  end {filter_mach};


procedure filter_psi;
  begin
    while psi > 2*PI do psi :=  psi - 2*PI;
    while psi < 0000 do psi :=  psi + 2*PI;
  end; { filter_psi }


procedure filter_data;
  begin
    filter_weight;
    filter_altitude;
    filter_rpm;
    filter_mach;
    filter_psi;
  end {filter_data};


procedure initialize_variables;
  begin
    GetAspectRatio(  Xasp,  Yasp );
    gAR   :=             Xasp/Yasp;

    col_1 :=         GetMaxX div 8;
    col_2 := col_1 + GetMaxX div 4;
    col_3 := col_2 + GetMaxX div 4;
    col_4 := col_3 + GetMaxX div 4;

    row_1 :=         GetMaxY div 6;
    row_2 := row_1 + GetMaxY div 3;
    row_3 := row_2 + GetMaxY div 3;

    radii := GetMaxX div 12; R := round( 0.9 * radii );

    ias_x := col_1;    hdi_x := col_2;    alt_x := col_3;    rpm_x := col_4;
    ias_y := row_1;    hdi_y := row_1;    alt_y := row_1;    rpm_y := row_1;
    ias_r := radii;    hdi_r := radii;    alt_r := radii;    rpm_r := radii;

    tab_x := col_1;    hsi_x := col_2;    roc_x := col_3;
    tab_y := row_2;    hsi_y := row_2;    roc_y := row_2;
    tab_r := radii;    hsi_r := radii;    roc_r := radii;

    attitude      :=          2.5 * PI/180;
    psi           :=        240.0 * PI/180;
    rpm           :=   75 / max_engine_rpm;
    mach          :=                 280.0;
    altitude      :=               25000.0;
    w             :=  0.9 * maximum_weight;
    theta         :=        0.001 * PI/180;
    alpha_tail    :=                   0.0;
    alpha_aileron :=                   0.0;
    trim          := m / (tail_arm * q * k_elevator);
    filter_data;
  end {initialize_variables};


procedure draw_rpm(x, y : integer);
  begin
    rpm_last  := rpm;
    str( rpm_last * max_engine_rpm:5:0, s );
    outtextxy(rpm_x - 5 * (text div 2), rpm_y + text div 2, s );
    line(rpm_x, rpm_y, x, y);
    rpm_lastx := x;
    rpm_lasty := y;
  end;


procedure rpm_needle;
  var   x, y   : integer;  z   : real;
  begin
    z := rpm * (max_engine_rpm/100) * 16 * PI/10;
    x := rpm_x + round( 0.9 * rpm_r * cos(13*PI/10 - z) );
    y := rpm_y - round( 0.9 * rpm_r * sin(13*PI/10 - z) * gAR );
    if  ( (abs(rpm_lastx - x) > 0)  OR  (abs(rpm_lasty - y) > 0) )
    then  begin
          setcolor(black);
          str( rpm_last * max_engine_rpm:5:0, s );
          outtextxy(rpm_x-5 * (text div 2), rpm_y + text div 2, s );
          line(rpm_x, rpm_y, rpm_lastx, rpm_lasty);
          setcolor(white);
          draw_rpm(x, y);
    end;
  end; { rpm_needle }


procedure rpmmeter;
  var z : real;
  begin
    circle(rpm_x, rpm_y, rpm_r);
    outtextxy(rpm_x+             round(rpm_r*sin((13-2*0)*PI/10)),
              rpm_y-             round(gAR*rpm_r*sin((13-2*0)*PI/10)),'0');
    outtextxy(rpm_x-(text  *  2)+round(rpm_r*cos((13-2*1)*PI/10)),
              rpm_y-             round(gAR*rpm_r*sin((13-2*1)*PI/10)),'12');
    outtextxy(rpm_x-(text  *  2)+round(rpm_r*cos((13-2*2)*PI/10)),
              rpm_y-(text div 2)-round(gAR*rpm_r*sin((13-2*2)*PI/10)),'25');
    outtextxy(rpm_x-(text  *  2)+round(rpm_r*cos((13-2*3)*PI/10)),
              rpm_y-(text div 2)-round(gAR*rpm_r*sin((13-2*3)*PI/10)),'37');
    outtextxy(rpm_x-(text      )+round(rpm_r*cos((13-2*4)*PI/10)),
              rpm_y-(text      )-round(gAR*rpm_r*sin((13-2*4)*PI/10)),'50');
    outtextxy(rpm_x+(text div 2)+round(rpm_r*cos((13-2*5)*PI/10)),
              rpm_y-(text div 2)-round(gAR*rpm_r*sin((13-2*5)*PI/10)),'62');
    outtextxy(rpm_x+(text div 2)+round(rpm_r*cos((13-2*6)*PI/10)),
              rpm_y-(text div 2)-round(gAR*rpm_r*sin((13-2*6)*PI/10)),'75');
    outtextxy(rpm_x+(text div 2)+round(rpm_r*cos((13-2*7)*PI/10)),
              rpm_y-(text div 2)-round(gAR*rpm_r*sin((13-2*7)*PI/10)),'87');
    outtextxy(rpm_x+(text div 2)+round(rpm_r*cos((13-2*8)*PI/10)),
              rpm_y-             round(gAR*rpm_r*sin((13-2*8)*PI/10)),'100');
  end; { rpmmeter }


procedure draw_n;
  begin
    n_last := n;
    str( n_last:5:2, s );
    outtextxy(ias_x - 6 * (text div 2),
              ias_y - round(gAR * ias_r) - text, s);
  end;  { draw_n }


procedure draw_mach;
  begin
    mach_last  := mach;
    str( mach_last:5:3, s );
    outtextxy( ias_x -     5     * (text div 2),
               ias_y - round(gAR * ias_r) + text, s );
  end;  { draw_mach }


procedure draw_tas;
  begin
    tas_last  := V;
    str( tas_last * feet_to_knots:3:0, s );
    outtextxy( ias_x - 3 * (text div 2 ),
               ias_y - round(gAR * ias_r) + 3 * text, s );
  end;  { draw_tas }


procedure draw_ias(x, y : integer);
  begin
    line(ias_x, ias_y, x, y);
    tas_lastx := x;
    tas_lasty := y;
  end;  { draw_ias }


procedure ias_needle;
  const
         range  =  500;
         k      =  5.5;  { ~3.5*PI/2  }
         e      =    0;
  var    x, y : integer;
  begin
    x := ias_x + round(0.8 * ias_r
            * cos( (ias * 5.5 / (range * knots_to_feet)) - PI/2) );
    y := ias_y + round(gAR * 0.8 * ias_r
            * sin( (ias * 5.5 / (range * knots_to_feet)) - PI/2));

    if abs(n - n_last) >= 0.01 then begin
      setcolor(black);
      str( n_last:5:2, s );
      outtextxy(ias_x - 6*(text div 2),
                ias_y - round(gAR * ias_r) - text, s);
      setcolor(white);
      draw_n;
    end;

    if abs(mach - mach_last) >= 0.001 then begin
      setcolor(black);
      str(mach_last:5:3, s);
      outtextxy( ias_x -     5     * (text div 2),
                 ias_y - round(gAR * ias_r) + text, s );
      setcolor(white);
      draw_mach;
    end;

    if abs(V - tas_last) * feet_to_knots >= 1 then begin
      setcolor(black);
      str( tas_last * feet_to_knots:3:0, s );
      outtextxy( ias_x -     3     * (text div 2 ),
                 ias_y - round(gAR * ias_r) + 3 * text,    s );
      setcolor(white);
      draw_tas;
    end;

    if  ( (abs(tas_lastx - x) > e)  OR  (abs(tas_lasty - y) > e) ) then
    begin
      setcolor(black);
      line(ias_x, ias_y, tas_lastx, tas_lasty);
      setcolor(white);
      draw_ias(x, y);
    end;
  end; {ias_needle}


procedure airspeed;
  const
    k_range = 0.011;  { arc(~3.5*PI/2) / max value 500  }
  begin {airspeed}              {  5.5 / 500  }
    circle(ias_x, ias_y, ias_r);
    outtextxy(ias_x + 3 * text div 2 + round(ias_r * cos((060*k_range) - PI/2)),
              ias_y + round(gAR*ias_r * sin((060*k_range) - PI/2)), '60');
    outtextxy(ias_x + 1 * text + round(ias_r * cos((120*k_range) - PI/2)),
              ias_y + round(gAR*ias_r * sin((120*k_range) - PI/2)), '120');
    outtextxy(ias_x + 1 * text + round(ias_r * cos((180*k_range) - PI/2)),
              ias_y + round(gAR*ias_r * sin((180*k_range) - PI/2)), '180');
    outtextxy(ias_x + 0 * text + round(ias_r*cos((240*k_range)-PI/2)),
              ias_y + round(gAR*ias_r * sin((240*k_range) - PI/2)),'240');
    outtextxy(ias_x - 2 * text + round(ias_r*cos((300*k_range)-PI/2)),
              ias_y + round(gAR*ias_r * sin((300*k_range) - PI/2)),'300');
    outtextxy(ias_x - 3 * text + round(ias_r*cos((360*k_range)-PI/2)),
              ias_y + round(gAR*ias_r * sin((360*k_range) - PI/2)),'360');
    outtextxy(ias_x - 3 * text + round(ias_r*cos((420*k_range)-PI/2)),
              ias_y + round(gAR*ias_r * sin((420*k_range) - PI/2)),'420');
    outtextxy(ias_x - 7 * text div 2 + round(ias_r*cos((480*k_range)-PI/2)),
              ias_y + round(gAR*ias_r * sin((480*k_range) - PI/2)),'480');
  end; { airspeed }


function Ydenom : real;
  var x : real;
  begin
    x := sqr( sin(attitude) * sqr(cos(theta)) );
    if abs(x) < dev then x := dev;
    Ydenom := x;
  end;


function Yw : integer;
  begin
    Yw := round( gAR * R * sin(attitude) * sqr( cos(theta) )
          * ( 1 + sqrt(
          abs( 1 + ( sqr(sin(theta)) - sqr(sin(attitude)*cos(theta)))
          / Ydenom ))));
  end;


function Yh : integer;
  begin
    Yh := round( gAR * R * sin(attitude) * sqr( cos(theta) )
          * ( 1 - sqrt(
          abs( 1 + ( sqr(sin(theta)) - sqr(sin(attitude)*cos(theta)))
          / Ydenom ))));
  end;


function Xdenom : real;
  var x : real;
  begin
    x := sqr( sin(theta) * sin(attitude) );
    if abs(x) < dev then x := dev;
    Xdenom := x;
  end;


function Xw : integer;
  begin
    Xw := round(-R * sin(attitude) * sin(theta) * cos(theta)
          * ( 1 + sqrt(abs( 1 + sqr( cos(attitude) ) / Xdenom ) ) ) );
  end;


function Xh : integer;
  begin
    Xh := round(-R * sin(attitude) * sin(theta) * cos(theta)
          * ( 1 - sqrt(abs( 1 + sqr( cos(attitude) ) / Xdenom ) ) ) );
  end;


procedure draw_hdi;
  var  x0, y0 : integer;
  begin
    sky_x := hdi_x - round( 0.8 * hdi_r * sin(theta) * cos(attitude)
                            / abs(cos(attitude)) );
    sky_y := hdi_y - round( 0.8 * gAR * hdi_r * cos(theta) * cos(attitude)
                            / abs(cos(attitude)) ) - text div 2;

    outtextxy(sky_x - (text div 2), sky_y, #30);

    xl := hdi_x + Xw;
    xm := hdi_x + Xh;
    yl := hdi_y + Yw;
    ym := hdi_y + Yh;


    line(xl, yl, xm, ym);
    line(hdi_x - hdi_r div 2, hdi_y + hdi_r div 6, hdi_x, hdi_y);
    line(hdi_x + hdi_r div 2, hdi_y + hdi_r div 6, hdi_x, hdi_y);
    line(hdi_x - hdi_r div 2, hdi_y + hdi_r div 6,
         hdi_x + hdi_r div 2, hdi_y + hdi_r div 6);
    circle( (xl+xm) div 2, (yl+ym) div 2, 4 );

    circle(hdi_x, hdi_y, 3);

    att_last   := attitude;
    theta_last :=    theta;
  end; { draw_hdi }


procedure horizon_bar;
  begin
    if (  ( abs(attitude - att_last) >= 0.005760 )   { about 0.33 deg }
       or ( abs(theta - theta_last)  >= 0.017453 ) ) { about 0.99 deg }
    then begin
         setcolor(black);
         outtextxy(sky_x - (text div 2), sky_y, #30);
         circle( (xl+xm) div 2, (yl+ym) div 2, 4 );
         line(xl, yl, xm, ym);
         setcolor(white);
         draw_hdi;
    end;
  end; { horizon_bar }


procedure hdi;
  var    s1, s2, s3, s4   : integer;
  begin
    s1 := round(hdi_r * sin( PI/6 ));
    s2 := round(gAR * hdi_r * cos( PI/6 ));
    s3 := round(hdi_r * sin( PI/3 ));
    s4 := round(gAR * hdi_r * cos( PI/3 ));
    circle(hdi_x, hdi_y, hdi_r);
    setlinestyle(solidln, 0, thickwidth);
    line(hdi_x-1, hdi_y-round(gAR * hdi_r),
         hdi_x-1, hdi_y-round(gAR * hdi_r)-text);
    line(hdi_x-s1, hdi_y-s2, hdi_x-s1-4, hdi_y-s2-4);
    line(hdi_x+s1, hdi_y-s2, hdi_x+s1+4, hdi_y-s2-4);
    line(hdi_x-s3, hdi_y-s4, hdi_x-s3-4, hdi_y-s4-4);
    line(hdi_x+s3, hdi_y-s4, hdi_x+s3+4, hdi_y-s4-4);
    setlinestyle(solidln, 0, normwidth);
  end; { hdi }


procedure draw_altitude(x, y : integer);
  begin
    alt_last  := altitude;
    str( alt_last:5:0, s );
    outtextxy(alt_x - 5 * text div 2, alt_y + text, s);
    line(alt_x, alt_y, x, y);
    alt_lastx := x;
    alt_lasty := y;
  end;


procedure alt_needle;
  const e      =       0;
  var   x, y   : integer;
        z      :    real;
  begin
    z := round( altitude / 1000 );
    z := round( altitude - (z * 1000) );
    x := alt_x + round( 0.82 * alt_r * cos(PI/2 - ((z/1000) * 2 * PI)));
    y := alt_y - round( 0.82 * alt_r * sin(PI/2 - ((z/1000) * 2 * PI))*gAR);
    if  ( (abs(alt_lastx - x) > 0)  OR  (abs(alt_lasty - y) > 0) )
    then  begin
          setcolor(black);
          str( alt_last:5:0, s );
          outtextxy(alt_x - 5 * text div 2, alt_y + text, s);
          line(alt_x, alt_y, alt_lastx, alt_lasty);
          setcolor(white);
          draw_altitude(x, y);
    end;
  end; { alt_needle }


procedure altimeter;
  begin
    circle(alt_x, alt_y, alt_r);
    outtextxy(alt_x - (text div 2),
              alt_y - (text      )-round(gAR*alt_r*sin((2*0+5)*PI/10)),'0');
    outtextxy(alt_x -              round(alt_r*cos((2*1+5)*PI/10)),
              alt_y - (text      )-round(gAR*alt_r*sin((2*1+5)*PI/10)),'1');
    outtextxy(alt_x + (text  -  2)-round(alt_r*cos((2*2+5)*PI/10)),
              alt_y - (text div 2)-round(gAR*alt_r*sin((2*2+5)*PI/10)),'2');
    outtextxy(alt_x + (text  -  2)-round(alt_r*cos((2*3+5)*PI/10)),
              alt_y -              round(gAR*alt_r*sin((2*3+5)*PI/10)),'3');
    outtextxy(alt_x -              round(alt_r*cos((2*4+5)*PI/10)),
              alt_y -              round(gAR*alt_r*sin((2*4+5)*PI/10)),'4');
    outtextxy(alt_x - (text div 2)-round(alt_r*cos((2*5+5)*PI/10)),
              alt_y + (text div 2)-round(gAR*alt_r*sin((2*5+5)*PI/10)),'5');
    outtextxy(alt_x - (text      )-round(alt_r*cos((2*6+5)*PI/10)),
              alt_y -              round(gAR*alt_r*sin((2*6+5)*PI/10)),'6');
    outtextxy(alt_x - (text  +  2)-round(alt_r*cos((2*7+5)*PI/10)),
              alt_y -              round(gAR*alt_r*sin((2*7+5)*PI/10)),'7');
    outtextxy(alt_x - (text  +  2)-round(alt_r*cos((2*8+5)*PI/10)),
              alt_y - (text div 2)-round(gAR*alt_r*sin((2*8+5)*PI/10)),'8');
    outtextxy(alt_x - (text      )-round(alt_r*cos((2*9+5)*PI/10)),
              alt_y - (text      )-round(gAR*alt_r*sin((2*9+5)*PI/10)),'9');
  end; { altimeter }


procedure draw_t_and_b(x, y : integer);
  begin
    line(tab_x, tab_y, tab_x + y, tab_y - x);
    setlinestyle(solidln, 0, normwidth);
    tbx := x;
    tby := y;
  end;  { draw_t_and_b }


procedure t_and_b_needle;
  const  k    =      10;
  var    t    :    real;
         x, y : integer;
  begin
    t := psidot * k;
    if abs(t) > PI/4   then    t := (PI/4) * t / abs(t);
    x := round( 0.83 * tab_r * cos(t) * gAR);
    y := round( 0.83 * tab_r * sin(t) );
    if ((tby <> y) or (tbx <> x)) then begin
      setlinestyle(solidln, 0, thickwidth);
      setcolor(black);
      line(tab_x, tab_y, tab_x + tby, tab_y - tbx);
      setcolor(white);
      draw_t_and_b(x, y);
    end;
  end; { t_and_b_needle }


procedure t_and_b;
  var    s1, s2 : integer;
  begin
    s1 := round(tab_r * sin( PI/6 ) );
    s2 := round(tab_r * cos( PI/6 ) * gAR);
    circle(tab_x, tab_y, tab_r);
    setlinestyle(solidln, 0, thickwidth);
    line(tab_x, tab_y-round(gAR*tab_r),
         tab_x, tab_y-round(gAR*tab_r)-text div 2);
    line(tab_x-s1, tab_y-s2, tab_x-s1-4, tab_y-s2-4);
    line(tab_x+s1, tab_y-s2, tab_x+s1+4, tab_y-s2-4);
    setlinestyle(solidln, 0, normwidth);
  end; {t_and_b}


procedure draw_heading;
  begin
    psi_last  := psi;
    str( psi_last * 180/PI:3:0, s );
    outtextxy(hsi_x - length(s) * (text div 2),
              hsi_y + text - round(gAR*hsi_r), s);
  end;  { draw_heading }


procedure hsi_heading;
begin
  if psi <> psi_last then begin
    setcolor(black);
    str( psi_last*180/PI:3:0, s );
    outtextxy(hsi_x - length(s) * (text div 2),
              hsi_y + text - round(gAR*hsi_r), s);
    setcolor(white);
    draw_heading;
  end;
end; { hsi_heading }


procedure hsi;
  begin
    circle(hsi_x, hsi_y, hsi_r);
    setlinestyle(solidln, 0, thickwidth);
    outtextxy(hsi_x-(text div 2), hsi_y-round(gAR*hsi_r), 'v');
    setlinestyle(solidln, 0, normwidth);
  end; { hsi }


procedure draw_rofc(x, y : integer);
  begin
    rcx := x; rcy := y;
    line( roc_x, roc_y, x, y);
  end;  { draw_rofc }


procedure rofc_needle;
  const  k     =  0.1132;    { 6000 ft/min = full scale i.e. 90%  }
  var    x, y  : integer;    t : real;
  begin
    t := hdot * k;
    if abs(t) > 0.9 * PI   then   t := 0.9 * PI * t / abs(t);
    x := roc_x - round(0.82 * roc_r * cos(t) );
    y := roc_y - round(0.82 * roc_r * sin(t)*gAR );
    if ( (abs(rcx - x) > 0) or (abs(rcy - y) > 0) ) then begin
      setcolor(black);
      line( roc_x, roc_y, rcx, rcy);
      setcolor(white);
      draw_rofc(x, y);
    end;
  end; { rofc_needle }


procedure rofc;
  begin
    circle(roc_x, roc_y, roc_r);
    outtextxy(roc_x - ( 3 * text div 2) - roc_r,
              roc_y-( text div 2),'0');
    outtextxy(roc_x-(text  *  3) - round(roc_r*cos( 0.9*1*PI/6)),
              roc_y-(text      ) - round(gAR*roc_r*sin( 0.9*1*PI/6)),'  1');
    outtextxy(roc_x-(text  *  3) - round(roc_r*cos(-0.9*1*PI/6)),
              roc_y-               round(gAR*roc_r*sin(-0.9*1*PI/6)),'  1');
    outtextxy(roc_x-(text      ) - round(roc_r*cos( 0.9*2*PI/6)),
              roc_y-(text      ) - round(gAR*roc_r*sin( 0.9*2*PI/6)),'2');
    outtextxy(roc_x-(text      ) - round(roc_r*cos(-0.9*2*PI/6)),
              roc_y-               round(gAR*roc_r*sin(-0.9*2*PI/6)),'2');
    outtextxy(roc_x-(text div 2) - round(roc_r*cos( 0.9*4*PI/6)),
              roc_y-(text      ) - round(gAR*roc_r*sin( 0.9*4*PI/6)),'4');
    outtextxy(roc_x-(text div 2) - round(roc_r*cos(-0.9*4*PI/6)),
              roc_y+(text div 2) - round(gAR*roc_r*sin(-0.9*4*PI/6)),'4');
    outtextxy(roc_x+(text      ) - round(roc_r*cos( 0.9*6*PI/6)),
              roc_y-(text      ) - round(gAR*roc_r*sin( 0.9*6*PI/6)),'6');
    outtextxy(roc_x+(text      ) - round(roc_r*cos(-0.9*6*PI/6)),
              roc_y-(text div 2) - round(gAR*roc_r*sin(-0.9*6*PI/6)),'6');
  end; {rofc}


procedure draw_fuelflow;
  begin
    last_fuelflow := fuelflow;
    str( (last_fuelflow * seconds_to_hour):4:1, s );
    outtextxy( col_4 - (6 * text div 2),
               row_1+(row_2-row_1-text) div 2, 'ff '+s );
  end;  { draw_fuelflow }


procedure fuelflow_value;
  begin
    if ( abs(last_fuelflow - fuelflow) > (1/6000) ) then begin
      setcolor(black);
      str( (last_fuelflow * seconds_to_hour):4:1, s );
      outtextxy( col_4-(6 * text div 2),
                 row_1+(row_2-row_1-text) div 2, 'ff '+s );
      setcolor(white);
      draw_fuelflow;
    end;
  end; { fuelflow_value }


procedure draw_circles;
  begin
    airspeed;
    draw_n;
    draw_mach;
    draw_tas;
    draw_ias(ias_x,ias_y);
    hdi;
    draw_hdi;
    altimeter;
    draw_altitude(alt_x,alt_y);
    t_and_b;
    draw_t_and_b(0,0);
    hsi;
    draw_heading;
    rofc;
    draw_rofc(roc_x,roc_y);
    draw_fuelflow;
    rpmmeter;
    draw_rpm(rpm_x,rpm_y);
  end; { draw_circles }


procedure print_data;
  begin {print_data}
    ias_needle;
    horizon_bar;
    alt_needle;
    t_and_b_needle;
    hsi_heading;
    rofc_needle;
    fuelflow_value;
    rpm_needle;
  end {print_data};


procedure calculate_time_changes;
  begin
    psi      := psi + psidot * dtime;
    theta    := theta + thetadot * dtime;
    attitude := attitude + gammadot * dtime;
    w        := w - dtime * fuelflow;
    altitude := altitude + hdot * dtime;
    mach     := calculate_mach_from_tas((V+Vdot * dtime)*feet_to_knots);
    filter_data;
  end {calculate_time_changes};


procedure readreal(var realnumber : real);
  const  bs = #8; { backspace }
  var    keys   : string;
         number : real;
         code   : integer;
         funckey: boolean;
  procedure decrease;
    begin
      if length(keys) > 0 then begin
        keys := copy(keys, 1, length(keys)-1);
      end;
    end; { decrease }
  begin
    keys    := '';
    repeat
      ch    := readkey;
      if ch <> #0 then funckey := false
      else
        begin
          funckey := true;
          ch      := readkey;
      end;
      case funckey of
        true  : case ord(ch) of
                  75  :  decrease;
                end;
        false : case ch of
                  '-','+','0'..'9', '.' : begin
                                            keys := keys + ch;
                                          end;
                end;
      end;
    until ch in [ #10, #13 ];
    if pos('.', keys) = 1 then keys := '0' + keys;
    val( keys, number, code );
    if  (code=0)  then   realnumber := number;
  end {readreal};


procedure set_rpm;
  var  x  : real;
  begin
    x := rpm * max_engine_rpm;
    readreal(x);
    x := x / max_engine_rpm;
    rpm  := x;
    filter_rpm;
  end; {set_rpm}


procedure set_pitch;
  var  x  : real;
  begin
    x := attitude * 180/PI;
    readreal(x);
    attitude := x * PI/180;
  end; {set_pitch}


procedure set_elevator;
  var  x  : real;
  begin
    x := alpha_tail * 180/PI;
    readreal(x);
    x := x * PI/180;
    if x >  max_alpha_tail then x :=  max_alpha_tail;
    if x < -max_alpha_tail then x := -max_alpha_tail;
    alpha_tail  :=   x;
  end; {set_elevator}


procedure set_aileron;
  var  x  : real;
  begin
    x := alpha_aileron * 180/PI;
    readreal(x);
    x := x * PI/180;
    if x >  max_aileron then x :=  max_aileron;
    if x < -max_aileron then x := -max_aileron;
    alpha_aileron := x;
  end; {set_aileron}


procedure set_mach;
  var  x  : real;
  begin
    x := mach;
    readreal(x);
    mach  :=  x;
    filter_mach;
  end; {set_mach}


procedure set_weight;
  var  x  : real;
  begin
    x := w;
    readreal(x);
    w := x;
    filter_weight;
  end; {set_weight}


procedure set_bank;
  var  x  : real;
  begin
    x := theta * 180/PI;
    readreal(x);
    theta  := x * PI/180;
  end; {set_bank}


procedure set_heading;
  var  x  : real;
  begin
    x := psi * 180/PI;
    readreal(x);
    psi :=  x * PI/180;
    filter_psi;
  end; {set_heading}


procedure set_altitude;
  var  x  : real;
  begin
    x := altitude/100;
    readreal(x);
    altitude := x * 100;
    filter_altitude;
  end; {set_altitude}


procedure set_level; { pitch is wrong in inverted flight }
  const  e = 0.09;   { ~cos(85) }
  begin
    if ( (abs(cos(theta)) > e) and (abs(cos(alpha)) > e) ) then
    begin
      alpha_tail := 0.0;
      trim       := m / ( Tail_arm * k_elevator * q * cos(theta) );
      rpm        := D / ( maximum_thrust * cos(alpha) );
      attitude   := alpha;
      filter_data;
    end;
  end; { set_level }


procedure set_trim;
  begin
    trim       := alpha_tail + trim;
    alpha_tail :=               0.0;
  end; { set_trim }


procedure idle;
  begin
    while not keypressed do
      begin
        calculate_time_changes;
        print_data;
      end;
    read_character;
  end; {idle}


procedure pickupkeys;  { keyboard codes on page 579 of Pascal manual }
  const delta_aileron  = 2.0E-3; { These are all }
        delta_elevator = 1.0E-3; { just a bunch  }
        delta_rpm      = 1.0E-2; { of guesses    }
  begin
    repeat
      idle;
      case function_key of
        true  : case ord(ch) of
        { PgUp }  73 : begin
                         rpm := rpm + delta_rpm;
                         filter_rpm;
                       end;
        { PgDn }  81 : begin
                         rpm := rpm - delta_rpm;
                         filter_rpm;
                       end;
        { UpArw } 72 : begin
                         alpha_tail := alpha_tail - delta_elevator;
                         if alpha_tail        <    -max_alpha_tail
                            then alpha_tail  :=    -max_alpha_tail;
                       end;
        { DwnArw }80 : begin
                         alpha_tail := alpha_tail + delta_elevator;
                         if alpha_tail         >    max_alpha_tail
                            then alpha_tail   :=    max_alpha_tail;
                       end;
        { LftArw }75 : begin
                         alpha_aileron := alpha_aileron - delta_aileron;
                          if alpha_aileron        <        -max_aileron
                             then alpha_aileron  :=        -max_aileron;
                       end;
        { RtArw } 77 : begin
                         alpha_aileron := alpha_aileron + delta_aileron;
                         if alpha_aileron         >         max_aileron
                            then alpha_aileron   :=         max_aileron;
                       end;
        { End }   79 : begin
                         alpha_aileron := 0;
                         alpha_tail    := 0;
                       end;
        { Home }  71 : begin
                         alpha_aileron := 0;
                         theta         := 0.001;
                         alpha_tail    := 0;
                         trim          := m/(tail_arm * q * k_elevator);
                         attitude      := 2.25 * PI/180;
                       end;
        {
         ^LArw = 115 : ^RArw = 116 : ^End  = 117 :
         ^PgDn = 118 : ^Home = 119 : ^PgUp = 132 :
        }
                end;
        false : case upcase(ch) of
                    'R' : set_rpm;
                    'P' : set_pitch;
                    'E' : set_elevator;
                    'A' : set_aileron;
                    'S' : set_mach;
                    'W' : set_weight;
                    'B' : set_bank;
                    'F' : set_altitude;
                    'H' : set_heading;
                    'L' : set_level;
                    'T' : set_trim;
                end;
      end;
    until ((ch in ['Q','q']) and (function_key=false));
  end; { pickupkeys }


{$F+}    (* force far call *)
procedure CleanUp;
  begin
    ExitProc := SaveExitProc;
    CloseGraph;
  end; { CleanUp }
{$F-}


procedure SetGraphics;
  begin
    { this may need modification if graphics do not conform to AT&T format }
(*****
    GraphDriver  := ATT400;     { use this for AT&T format }
    GraphMode    := ATT400Hi;   { 2 color 640 X 400  }
*****)

    GraphDriver  := detect;     { use this to autodetect graphics format }

    InitGraph( GraphDriver, GraphMode, '');
    ErrorCode    := GraphResult;
    if ErrorCode <> grOk then begin
      Writeln('Graphics error: ', GraphErrorMsg(ErrorCode));
      writeln('Program aborted...');
      halt(1);
    end;
  end; { SetGraphics }


Begin {Simulate}
  clrscr;
  SaveExitProc := ExitProc;
  ExitProc     := @CleanUp;
  SetGraphics;
  setviewport(0, 0, getmaxx, getmaxy, clipon);
  initialize_variables;
  set_level;

  repeat
    clearviewport;
    draw_circles;
    print_data;
    pickupkeys;
    moveto( (GetMaxX - 27 * 8) div 2, GetMaxY - 8 );
    OutText('Quit = Q : <cr> = Continue ' );
    read_character;
  until ch in ['q','Q'];

End {Simulate}.


      (***********************************************************
      This software is made available only to individuals and only
      for educational purposes. Any and all commercial use is
      stricly prohibited.
      ***********************************************************)




      (************************   END   **************************
      **                                                        **
      **                                                        **
      **                    Flight  Simulator                   **
      **                  Subsonic Jet Aircraft                 **
      **                      Version  4.26                     **
      **                                                        **
      **                                                        **
      ************************************************************

      ************************************************************
      *                                                          *
      *                                                          *
      *            Captain Larry Thomas Brewster, Ph.D.          *
      *                                                          *
      *                University of Missouri-Rolla              *
      *                Computer Science Department               *
      *                       Feb 26, 1990                       *
      *                                                          *
      *                                                          *
      ***********************************************************)